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ABSTRACT 

Anomalous X-ray Pulsars (AXPs) and Soft Gamma-Ray Repeaters (SGRs) are young neutron stars (NSs) 
characterized by high X-ray quiescent luminosities, outbursts, and, in the case of SGRs, sporadic giant flares. 
They are believed to be powered by ultra-strong magnetic fields (hence dubbed magnetars). The diversity 
of their observed behaviours is however not understood, and made even more puzzling by the discovery of 
magnetar-like bursts from "low-field" pulsars. Here we perform long-term 2D simulations that follow the evo- 
lution of magnetic stresses in the crust; these, together with recent calculations of the breaking stress of the 
neutron star crust, allow us to establish when starquakes occur. For the first time, we provide a quantitative 
estimate of the burst energetics, event rate, and location on the neutron star surface, which bear a direct rel- 
evance for the interpretation of the overall magnetar phenomenology. Typically, an "SGR-like" object tends 
to be more active than an "AXP-like" object or a "high-Z? radio pulsar", but there is no fundamental separa- 
tion among what constitutes the apparent different classes. Among the key elements that create the variety of 
observed phenomena, age is more important than a small variation in magnetic field strength. We find that 
outbursts can also be produced in old, lower-field pulsars (B ~ a few x 10 12 G), but those events are much less 
frequent than in young, high-field magnetars. 

Subject headings: stars: neutron — X-rays: stars 



1. INTRODUCTION 

Among the isolated NSs, particularly distinctive is a class 
of sources characterized by long periods (P ~ 2 — 1 1 s) and 
high quiescent X-ray luminosities (L x ~ 10 33 - 10 35 erg s -1 ), 
generally larger than their entire reservoir of rotational en- 
ergy. These sources, historically classified as AXPs and 
SGRs, often display stochastic bursts of X-rays, releasing en- 
ergies ~ 10 39 - 10 41 erg, and, in the case of SGRs, sporadic 
though very energetic 7-ray flares, with typical energetics 
~ 10 44 -10 45 erg. 

The most successful model to explain both the high qui- 
escent X-ray luminosities, as well as the X-ray bursts and gi- 
ant 7-ray fla res, is the magnetar model dThompson & Duncan! 
1 19951 120011) . in which SGRs and AXPs are believed to be en- 
dowed with large magnetic fields, B ~ 10 14 - 10 15 G, possi- 
bly resulting from an active dynamo at birth. After a mag- 
netar is born, the internal magnetic field is subject to a con- 
tinuous evolution through the processes of Ohmic dissipation, 
ambipolar diffusion, and Hall drift. In the crust, magnetic 
stresses are generally balanced by elastic stresses. However, 
as the internal field evolves, local magnetic stresses can oc- 
casionally become too strong to be balanced by the elastic 
strength of the crust, which hence breaks, and the extra stored 
magnetic/elastic energy becomes available for powering the 
bursts and flares. 

Despite the success of the magnetar model in explaining 
some general features of the triggering mechanism of bursts 
and flares, some major questions have been left unanswered. 
In particular: what determines the frequency of the bursts? 
Why do some objects display giant flares (SGRs), while oth- 
ers (AXPs) do not? How frequent are the different phenom- 
ena? Why do burst locations on the NS surface generally cor- 
relate with the pulsar phase (maximum of the quiescent X-ray 
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lightcurve) for AXPs, while they do not for SGRs? (see e.g. 
Kaspi & Gavriill (120041) for a review of these properties). 

Moreover, the discovery of magn etar-like X-ray burst s from 
the young pulsar PSR J 1846-0258 dGavriil et alj [2008). with 
an inferred surface dipolar magnetic field B p = 4.9 x 10 13 G, 
lower than the traditionally considered magnetar range, and, 
more recently, the discovery of SGR 0418+5 729 with an 
even lower B p < 7.5 x 10 12 G dRea et al.ll20Toh . well within 
the range of the rotation powered pulsars which do not dis- 
play any bursting behaviour, has raised another obvious ques- 
tion: why some "high-B" pulsars (PSR Jl 119-6127 and PSR 
J1814-1744) do not display any discernable X-ray emission 
nor outburst (C amilo et alJl2000l) . while at least one case of 
"low-B" NS does, if the magnetic field is their driving force ? 
It is clear that the dipolar magnetic field alone is not sufficient 
to account for this variety of behaviours. A unified physical 
framework is still lac king, and it c onstitutes a major puzzle in 
neutron star physics (Kaspi 2010). 

Here we present the first investigation that combines long- 
term 2D simulations of the coupled magneto-thermal evolu- 
tion of the NS with the study of the breaking of the NS crust. 
Our study, by using realistic values of the breaking stress of 
the NS crust, allows us to estimate burst energetics, recur- 
rence times, and surface distribution. Our approach allows us 
to identify some key elements of the magnetar phenomenol- 
ogy, and to shed light on what creates the variety of the - often 
puzzling - observed behaviour. 

2. MODELING THE COUPLED MAGNETIC-THERMAL EVOLUTION 
OF NSS AND CRUSTAL FRACTURES 

We follow the evolution (in axial symmetry) of the mag- 
netic field in a magnetar crust w ith the numerical code de- 
scribed in lPons & Geppertl (120071) . As in previous works, we 
restrict our study to magnetic field configurations confined to 
the crust. We note that one of the main products of our calcu- 
lations, i.e. heating by current dissipation, has little relevance 
in the core where the conductivity is high and the heat re- 
leased is lost by neutrino emission. Also important for our 
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work, the changes in elastic stresses and their effect in the 
crust are of no relevance in a liquid core. The only issue that 
may affect our results is the global change of magnetic field 
geometry that a fully coupled crust-core evolution may give, 
but the way in which this could change our picture is not clear. 
This shortcoming is caused by the difficulties (both theoreti- 
cal and numerical) associated to the superconducting state of 
protons in a NS core, which makes it difficult to consider the 
full problem at present, although we hope to be able to do so 
in future work. 

The initial configurations of the magnetic field include both 
a toroidal and a poloidal component. The temper ature of the 
crust a t different ages is given by the results of iPons et alJ 
(2009), where the interplay between the thermal and mag- 
netic field evolution of magnetars was studied. For simplic- 
ity, and for numerical limitations, here we assume that the 
crust is isothermal. As the magnetic field evolves, the crust 
moves through a series of equilibrium states in which its elas- 
tic stress <Tb(r,9,t) balances the (time-dependent) magnetic 
stress Mij(r,6,t) in each direction. Assuming that equilibrium 
is reached at a certain time, t eq , we can define 



4ti 



>a b (r 7 6,n, (1) 



where r and indicate the radial and poloidal coordinates, 
respectively. 

Recent molecul ar dynamical simulations 

(Chug unov & Horowitz! 1201 Oh have provided a fit for 
the maximum stress that a NS crust can sustain 
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where T = Z 2 e 2 /aT is the Coulomb coupling parameter, a = 
[3/(47m,-)] 1 ' 3 is the ion sphere radius, n,- the ion number den- 
sity, Z the charge number, e the electron charge, and T the 
temperature. 

At some time t b during the evolution, any component ij of 
the local magnetic stresses can occasionally depart from the 
previous equilibrium condition at f eq by an amount which is 
comparable to, or exceeds, the breaking stress of the crust, 
cr™ ax (r,fi) (angle-independent because of our assumption of 
isothermality); hence the crust fractures, yielding a starquake, 
and a new equilibrium state is reestablished. 

Since our simulations follow the evolution of B(r,0,t), and 
hence Mjj(r,0,t), we can map the time-dependent regions in 
the NS crust which are subject to fractures. A fully consistent 
dynamical simulation of the starquake is out of our present 
capabilities. The typical timescales (mseconds to seconds) 
of individual burst/flare events are many orders of magni- 
tude smaller than the long-term evolution timescales (typical 
timesteps in our code are ~ a week, to allow us to follow 
the NS evolution up to 10 5 — 10 6 years), hence we cannot dy- 
namically follow the fracture propagation or model individual 
bursts. However, we can estimate the energy of an "outburst" 
(i.e. a collection of tens to hundreds of bursts occurring within 
~ a week timescale) as follows. When a starquake happens, 
say at time t = ?/,, a certain region of the crust in the vicinity 
of the point where critical conditions are reached will be af- 
fected. We consider that all surrounding regions where local 
magnetic stresses are close to the maximum (say by a fraction 
e) are affected, i.e. 

M u (r, 0, t b )-M^(r, 0) > ea^(r, t„) . (3) 

This parameter e must be close to the fatigue limit of the 
material (for terrestrial materials it ranges between 35% and 



60%). This limit is also close to the yield point: the stress at 
which a material begins to deform p lastically. All investigated 
crystallographic shear systems in iHorowitz & Ka dau (2009) 
break in a rather abrupt fashion with only a small region where 
plasticity is present. Thus, we expect the parameter e to be in 
the range [0.8-0.99]. It parametrizes our ignorance about the 
transition from elastic to plastic regimes and the fatigue of the 
material. 

We can compute the elastic energy stored in that portion 
of the crust, and assume that it will be released at once (one 
timestep) and that the affected region will return immediately 
to equilibrium. The new equilibrium stresses are reset, and 
the process is repeated. We neglect the feedback produced 
by the local deposition of energy, which will require a much 
more complex modeling. In this way, we estimate the energy 
available in a given event, the time interval between events, 
and the location of the fractures. It should be stressed that 
this represents the total energy of an "outburst", that may be 
released in one very energetic flare, many small bursts, or a 
combination of both. Also, depending on the local physical 
conditions, part of the energy is lost to neutrinos, part is trans- 
ferred to the magnetosphere, and part results in local heating 
and is radiated by photons from the surface on a much longer 
timescale. Typically, energy released in the inner crust is more 
easily lost to neutrinos while energy released near the surface 
has a more direct observational impact through thermal (sur- 
face) and non-thermal (magnetospheric) radiation. 

The energy available after a starquake can therefore 
be estimated as the elastic energy corresponding to a 
static shear strain 'F, which can be well approximated by 
( Tho mpson & Duncanll2001l) 



E b (t b ) = 2 / dVfiw = E y / dV 



(4) 

where we have used the fact that the yield strain is approxi- 
mated by 



Mij(.r,e,t b )-M%(r,0) 



(5) 



and t he shear modulus /i is given by (IHorowitz & H ughto 
2008) 
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The integral in Eq. dU is computed over the volume dV for 
which condition (01 on the stresses in any direction is satis- 
fied. Additional energy is stored in the magnetic field, but the 
fraction of this energy released during a fracture is unknown, 
and depends on the amount of slippage. This requires a proper 
calculation involving fracture dynamics on short timescales. 
To be conservative, we consider the case in which only elastic 
energy is released, so our energy estimates are lower limits. 

3. PREDICTING STARQUAKE ENERGETICS, FREQUENCIES AND 
DISTRIBUTION ON NS SURFACE 

We now report the most relevant results. As a baseline 
model we have chosen a NS born with an initial dipolar 
field of B p = 8 x 10 14 G and an internal toroidal field of 
B, = 2 x 10 15 G at maximum, similar to models used in pre- 
vious studies and close to what one expects to be the final 
geometry a fter MHD equilibrium is reached in a hot, liq uid, 
proto-NS (ICiolfi et al.ll2009l l2010t Lander & Jonesll2009l) . In 
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this particular model, at the age of 10 4 yrs the dipolar field has 
decreased to about 20% of its initial value, while its internal 
toroidal field has been dissipated by more than a factor of 2. 
Other initial conditions modify quantitatively the results, but 
the general trends remain nearly the same. The evolution is 
followed for 10 5 years. During this time, we monitor the fre- 
quency, angular distribution, and energetics of the fractures. 




Fig. 1. — Snapshot of the deviation of stresses from the previous 
equilibrium in a typical 2000 yr old magnetar crust. The top, middle 
and bottom panels show the (8(j>), (r0), and (r8) components of the 
tensor My -M° q normalized to the local value of er™ ax . The crustal 
region has been stretched by a factor of 6 for clarity of visualization. 
The components of the stress tensor are normalized to the local value 
of the breaking stress, so that yellow regions are close to the fracture 
limit. 

Fig. 1 shows a snapshot of the deviation of stresses from the 



previous equilibrium in a typical 2000 yr old magnetar. The 
yellow regions are close to the breaking limit. Note that at 
this age the crust has already suffered hundreds of fractures. 
That history results in the patched appearance. For this run 
we have fixed e = 0.9. 

The results about frequency and energy distribution of 
events are contained in Fig. 2. We remind again that groups 
of bursts connected to the same fracture are classified as a 
single outburst in our simulations. We selected three repre- 
sentative periods in a magnetar life, each involving the same 
total number of events (1000). The first is labeled "Young" 
and spans the interval between 400 and 1600 yrs. The second 
is labeled "Mid age" and covers the period « 7 — 10 kyrs, and 
the third period is named "Old" and corresponds to the events 
recorded from 60 to 100 kyrs. Note that this classification 
does not pretend to reflect an exact correspondence to SGRs 
and AXPs, which is just terminology based on historical rea- 
sons. Objects that could belong to the first category (young or 
"SGR-like") are SGR 1806-20, SGR 1900+14 and IE 1547.0- 
5408; to the Mid Age or "AXP-like" list could belong SGR 
1627-41, SGR 0526-66, SGR 0501+4516, 1E1048.1-5937, or 
IE 2259+586 among others; finally, examples of old objects 
could be SGR 0418+5729, SGR 1833-0832, XTE J1810-197, 
and PSR J1814-1744. 

The first important result is that there is a significant differ- 
ence in the energetics and recurrence time as the star evolves. 
This is because the initial evolution of the crustal magnetic 
field is faster due to the Hall drift, while at late times its has 
been rearranged into a more quasi-steady state and the evolu- 
tion proceeds slower. In young magnetars, the typical ener- 
gies released after each starquake are of the order of 10 44 erg 
and the typical event rate ~ 1/yr. This does not necessarily 
imply that a flare will be observed, since this is the energy re- 
leased in the crust and the mechanism to transfer energy to the 
surface and magnetosphere may not work in many cases, es- 
pecially when the fracture occurs in the inner crust. However, 
in our simulations we observed that more than 95% of the 
fractures occur in the outer crust, simply because the crust is 
less strong there, so that we expect a relatively high efficiency 
of the process. For mid-age objects the energetics is shifted 
to lower values, a second peak appears at about 10 41 erg, and 
the waiting time between outbursts increases to a few years. 
This trend is more pronounced as the object gets older: the re- 
currence time becomes of the order of tens of years and nearly 
half the events are low-energy events. The physical reason be- 
hind the bimodal distribution is the directionality of stresses. 
We found that fractures associated to the magnetic stress com- 
ponent Me ,0 are more frequent, but they are mostly associated 
to very low energy events < 10 41 erg. On the other hand, frac- 
tures caused by the M r> ^ component are responsible for most 
of the E > 10 44 erg events. The events associated to large 
values of M r g span the whole range of energies, but they be- 
come very rare after a few kyrs, so that the long-term energy 
distribution becomes more clearly bimodal. On average, we 
find that events related to stresses created by the toroidal field 
are a hundred times more frequent. Our simulations hence 
predict that giant flares are expected to be less frequent than 
the less energetic bursts. Indeed, while giant flares have been 
observed in only 3 objects (once each), bursts have been ob- 
served in more than a dozen (often multiple times each). 

In Fig. 3 we show the results for the same case, but with the 
parameter e varying randomly in the interval [0.95-1] for each 
event. There is no reason why this parameter must be constant 
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and in fact fracture dynamics cannot be parametrized in a sim- 
ple form. Nevertheless, with this other study, we can assess 
the sensitivity of our results to the particular choice of e. The 
first important remark is that the energy distribution changes 
and now low-energy events are more frequent than giant-flare- 
like events even for young magnetars. This choice of values 
of e closer to unity reduces the average size of the region im- 
plied in the starquake, and selects less energetic events. Inter- 
estingly, the bimodal shape of the energy distribution remains 
qualitatively the same. In the central panel we can see that the 
age trend discussed in the previous case does not change: the 
older the object, the longer the average waiting times between 
two events. The comparison between these two cases indi- 
cates that a detailed understanding and modeling of fracture 
mechanics will be needed before we can make definite, quan- 
titative statements about the energy distribution of starquakes. 

It is interesting to compare the general trends predicted 
by our simulations with the observations. We find that, for 
a given fi-field strength and configuration at birth, younger 
objects are generally more active, both in event frequency 
and energy of events. Interestingly, this age sequence is als o 
what observations have been hinting at (lKouveliotoulll999h . 
Our predicted clustering of the energy distribution (around 
10 40 - 10 41 and 10 44 erg) is suggestive of the observed di- 
chotomy outbursts/giant flares, with rarer events in betwee n 
(iKouveliotou et al.ll200 lh llsrael et all2008tlMereghettil2008h . 
Note that our results predict that middle age objects also show 
high energy events; interestingly, some o f the originally clas- 
sified AXPs also show SGR-like bursts dKumar & Safi-harbl 
2010!). A rigorous comparison with observations cannot be 
made at this stage, due to the very limited statistics. However, 
if the secular waiting time distribution of starquakes (which 
is what we predict) has a similar behaviour to the waiting 
time distribution of bursts within an outbu rst, then our re- 
sults would agai n match the observations (Cheng Tt alJI 19961 
iGogus et aUlioOlh . 

Further connections to the observations can be made by in- 
spection of the right panels of Figs. 2 and 3, where we show 
the angular distribution of the bursts. While for "young" and 
"old" magnetars there is no clear trend, it is interesting to note 
that at intermediate age most of the fractures happen at small 
polar angles. If the polar region is hotter than the equator, our 
results suggest a correlation between burst location and pul- 
sar ph ase, which is typical of objects routinely classified as 
AXPs (iKaspi & Gavriill2004lG~avriil et alj2004b . The reason 
for this particular feature is probably related to the evolution 
of the internal toroidal field. The Hall drift displaces the inter- 
nal toroidal field towards the poles in a typical Hall timescale 
of 10 3 - 10 4 yrs, but after a few Hall timescales the magnetic 
field is reconfigured in a more stable, steady-state. The char- 



acteristic location of starquakes closer to the pole seems to be 
reflecting this fact. To confirm this point, consistent simula- 
tions that include the local heating and the temperature varia- 
tions produced by the energy release are in progress and will 
be presented elsewhere. 

Finally, we also studied the magnetic evolution of a NS with 
a lower dipolar magnetic field strength B p = 2 x 10 14 G, and a 
toroidal component B, = 10 15 G (at birth). We made a longer 
simulation covering 10 6 yr of this NS and found that, while 
the energetics is very similar, the event rate is much lower. 
Even at very early times the typical lapse time between events 
is about 50-100 yrs. At late times, it becomes larger, ranging 
from few kyrs to 10 kyrs. The reason why the energetics is 
similar is simple: the local physical conditions that establish 
the breaking criteria are the same (Eq. [3}, so the elastic en- 
ergy stored is similar. However, the magnetic field evolves 
slowly, so it takes more time to bring the system close to the 
critical conditions. Interestingly, our magneto-thermal evo- 
lution model naturally predicts that all NSs, not only magne- 
tars, may show activity in the form of sporadic outbursts or 
flares. This last case, for which B p becomes < 7 x 10 12 G at 
an age > 1 Myr, would be a cl ose representation of the burst- 
ing source recently reported bv lRea et al.l d2010h . Other cases 
studied show similar results, with the particularity that objects 
with a lower toroidal field are generally less active. 

4. SUMMARY 

In our study, we have found that there is no real separa- 
tion among what constitutes the apparently different classes of 
bursting NSs, and we have identified some key elements that 
create the variety of observed burst/flare phenomenology. In 
particular, while what we infer from measurements of P and 
P is only the dipolar component of the fi-field, both the dipo- 
lar component and the hidden toroidal component are simi- 
larly important. The lower the field (in either component), the 
lower is the frequency of starquakes. Although very rare, we 
find that outbursts can also occur in "low-B" pulsars, with B ~ 
a few x 10 12 G. For a given B field at birth, our results show 
that the age (or the evolutionary stage) is the most relevant 
feature to determine the different levels of activity in different 
families of NSs. The same object may behave differently at 
different times, without in principle having a different nature. 
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MEC grant AYA 2010-21097-C03-02, GVPROMETE02009- 
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the ESF (JP). We thank S. Mereghetti, N. Rea and an anony- 
mous referee for useful comments on the manuscript. 
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Fig. 2. — Outburst properties of an object with initial magnetic field components B p = 8 x 10 14 G and B, = 2 x 10 15 G during three different 
periods of its lifetime: 400-1600 yr (labeled "Young"), 7-10 kyr (labeled as "Mid age"), and 60-100 kry (labeled as "Old"). We have fixed 
e = 0.9. 



Energy Distribution 



3qq 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 300 



o 
o 
o 

- 200 



in 

£ 100 
o 



J1 



i ijjr^fr 



! i 
i 



o 

39 40 41 42 43 44 45 
log (E (erg)) 



Waiting time Distribution 



200 



100 







I I I I | I I I I | I I I I | I I I I | I I I I | I I I I 



Young 



Mid age 



Old 



u -' 

I 

_ I 

I ■ I ^.r-i i 



-1.0 -0.5 0.0 
log( 



i 



5 1.0 1.5 2.0 

,<y0) 



Starquake Location Distribution 
500 F -1 — ' — ' — ' — 1 — ' — ' — ' — ' — 1 — ' — ' — ' — ' — r~ 



400 



300 



200 



100 




0.5 1.0 
Angle (rod) 



Fig. 3. — Same as Fig. 2, but with e varying randomly in the interval [0.95-1]. 



